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ABSTRACT 


The earth tidal tilt data consists of two parts: one is 
the theoretical part caused by the tidal potential of the 
moon and the sun; another one is the seismic noise part that 
originates from meteorological, hydrodynamic, and cultural 
sources. The theoretical part is completely determinable 
by astronomical data. The seismic noise may be considered 
as a stationary stochastic process over short periods of 
time when there have been no great changes in the source 
mechanisms. Based on the above assumption, prediction and 
interpolation of the earth tidal tilt data can be achieved 
by solving the Wiener-Hopf equation. A computer program 
"PREDICT" is derived by using Levinson's recursive method to 
solve the discrete-time Wiener-Hopf equation. Test results 
on real data are satisfactory. Besides the earth tidal tilt 
data, the "PREDICT" program can also be applied to any 


stationary stochastic discrete-time signal. 
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I. INTRODUCTION 


The essential problem of earthquake prediction is to 
identify and monitor those physical characteristics of the 
failure process that uniquely herald the magnitude, location 
and origin time of an earthquake. Identification of these 
Saaracteristics has tended to be heuristic and statistical 
imoeapproach; however, evidence is steadily accumulating 
toward the support of a sound physical basis to earthquake 
prediction. Models of this failure process are becoming 
Suificiently detailed to effectively define the spatial and 
temporal character in each of the variety of geophysical 
measurements conducted in the epicentral region, 

This paper contains the theory and applications of a 
Single-channel predictor for tilt at tidal sensitivity. 
Tidal tilt is caused by the direct tidal forces of moon and 
Sun on the earth. There are also other sources such as 
meteorological activity and the ocean-loading effect which 
can cause significant tilts exceeding the magnitude of the 
primary source. The primary sources are predictable and, 
to a lesser degree, the effects of ocean and meteorological 
Memivity are also predictable. OGthexr factors such as effects 
local Carthquake activity are random in nature and cannot 
be determined in advance. Therefcre, we divide tidal-tilt 
Piornacton Into two classes: one is predictable, the other 


random, In the nomenclature of communication theory, we 
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call the former signals and the latter noises. However in 
Pwocequent discussion we will, in fact, use the definition 
that the signal we seek is random in nature and the noise 

is the periodic but undesired part of the time series, 

A single-channel predictor in time series analysis uses 
fest data to predict future data. The formulation of a self- 
predictor is basically one of deSigning a filter whose output 
in the time domain is a positive time-shifted version of the 
mepuc, the prediction error is a measure of the stationarity 
of the time series. If the time series spectral density of 
equal lengths of data segments is a constant, then a filter 
designed from any segment of past data can also be used to 
predict future data from any other data segment. For the 
geeophysical problem at hand, we derive such a filter during 
earthquake-free periods, because we assume that the failure 
process is, mathematically speaking, associated with non- 
mratwonarity in the data. The strategy is thus one of 
femavying a filter free of this "signal"; using that filter 
for subsequent incoming data, and monitoring the magnitude 
oa Che prediction error aS a means of identifying onset of 


the failure process. 





II, THEORY OF SINGLE-CHANNEL PREDICTOR 


A. BASIC TERMINOLOGIES 

mee Dicital Signals 

In order to use high-speed digital Pee to 

analyze seismic signals, the first thing we must do is to 
convert the analog signals to digital signals. This involves 
the work to extract the values of the analog signal at 
discrete time points which are usually equally svaced. It 
is convenient to choose the time unit equal to the time- 
spacing At between two consective discrete time points so 
that the discrete time variable t takes on only integer 
values. The numerical values of the signal together with 
their associated discrete time points make up a digital 


Signal. One example of a digital signal is 


puomiiemecre ft, 3) 2. 4°53, 5, 2, 4, +08 


iene eee a Oo. 1, 2 


? 


Se eck 


’ 


Because we have required time t to assume consecutive integer 


values, we can write the above digital signal simply as 


where the arrow is used to denote the time origin t = 0 and 


x represents the entire signal. Under this notation, xy 
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Mimeecdesignate the amplitude of the digital signal x at 


eome tt. 


2. 


Thus in the above signal Xo = 2, X, = 4, and x -— 3. 


Wavelets 


Wavelets are a special class of digital signals. A 


wavelet b, is characterized by two properties, namely, 


fe 
(a) The one-sided property: A wavelet has a definite 
Orleinmeer Arrival) time in the sense that all values 
Gieiemen wavelet before 1tS Origin» time are zero, that 
is by =e LOLet < 0’, 
(b) The stability property: A wavelet has finite energy 
so that it 1s a transient or dying-out phenomenon, 
Pat S 
ie 2 
> b, 
t=0 
is finite. 
From the above two properties we know that a wavelet 
is a one-sided energy Signal, It may damp out completely, 


1.e., become zero, after a certain time. This is the so- 


called 


aa 


finite-length wavelet. One example is 


, 0, 0, 0, 5, 2, 3, 0, 0, 0, “"*) 


where we can Simply write as 


(9, 


LIOR 


We call it a three-length wavelet. <A wavelet can also be of 


infinite-length, that is, it damps out at infinite time. One 


example of infinite-length wavelet is 


(1, 


o. 


ere eee cee loo t/ 32. 1/64, 1/128, 1/256, *""). 
Convolution in Discrete-Time Systems 


inmonderertOulind a linear System's response or output 


resulting from a specific input signal, we have two methods 


el 





in general use. One method is to model or characterize 
the linear system by means of linear differential 
(continuous-time systems) or difference (discrete-time 
systems) equations. The output of the system to any input 
Signal is then obtained by solving the appropriate equation. 
Another method is based on the superposition property of 
linear systems, using the convolution operation between the 
input signal and the system's impulse response to get the 
output signal. Since we are only interested in treating 
Seenals in the discrete-time domain, we consider the convo- 
lution operation here only in discrete-time systems. 

The convolution of two wavelets a, and b, is denoted 


1G 
by 


To illustrate the use of this formula let us suppose 


99 
| 


(ag, Qs: 25) Clive 2) and 


b = (b b b 


Ome a Nes 9) (25 1 A). 


Using the convolution formula, we have 


= K = 
c ) b (Co; Cy» Cor Cas Cy) 


G@, =a2a,b, ta_b. = 1l1l_xits3sx2e= 7 
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ae pete woe ano. =~ | x 44+43xi42+2x2= 11 


Z O°2 5 lial | 2 O 
Cg = abo + aob, =3x4+2xi1e= 14 
Cy = A5bo = 2x 4= 8 
Hence 


e 
I 

& 
% 

oO 
I 


Opera deed 8). 


In general case, if a, and b, are not wavelet but any two 


v 


digital signals, the formula for convolution operation is 


cHea* b 
where 
(so @) 
Cop ee een 
m=—0 


4. The Z-Transform 

In analyzing linear systems, we can use linear 
differential or difference equations. We can also use 
convolution operations. Both are called time-domain tech- 
niques since the independent parameter is time t. But we 
find that if we transform the time signals within our system 
to another form, we can often express relationships more 
Simply than we have previously in the time-domain. We call 
this transform approach transform-domain techniques. One 
important consequence of the transform-domain approach of 
linear systems is that the convolution operation in the 
time-domain is converted to a multiplication operation in 
the transform-domain. There are three transforms which are 
important in the study of linear systems. The Laplace 


transform and Fourier transform are used in the continuous-time 
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Syscems,. The Z-transform is appropriate for discrete- 

time systems. We will only consider the Z-transform here. 
The Z-transform of finite wavelet is the polynomial in Z 
whose coefficients are the coefficients of the wavelet. For 


example, the Z-transform of the wavelet 


a= (a5; Ay, Anyees, a) 
is the polynomial 
m4) = a.+ a.Z4ta z7 te Ao 
O 1 2 aa dl 


The Z-transform of an infinite wavelet is the power 
series in GZ whose coefficients are the coefficients of the 


wavelet. For example, the Z-transform of the wavelet 


b= (bo> bi» be ...-) is the power series 


7 2 
3(Z) S15 2: |e hy Ae a 


The Z-transform of a general two-sided digital signal 
which is of infinite length in both directions will be a 
Laurent series in Z. For example, the Z-transform of the 


digital signal 


C= CC... Cig» ©_4> fq» Cy) Co - ) 
is the Laurent series 
_ -2 -1 2 
C(Z) =... + Bf 1 OH eee c 2 + Co LZ. amen: 


To show convolution operation in time-domain may be 
performed by multiplication operation in transform domain, 


let us redo our example in the last section. We have 
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a= (ag, ay: ao) =e 3, 2) and 
a (bp, bi, bo) eee) 4 ) 
The Z-—-transforms of ay and bi are 
2 
A(Z) = 1+ 32+ 2424 and 
Paaye= 2 + Ze 47-- 


Multiplying A(Z) and B(Z), we obtain 


C(Z) = A(Z)B(Z) = 2 + 7Z + 112° + 14Z° + 82%. 


After we transform back to time-domain, we get 
c= (Co, Cy, Co, Ca, Cy) eZee iw ld 8) 


which is the same result as in the last section. 
Oo. Autocorrelation and Cross-Correlation 
Autocorrelation and cross-correlation are two 
important concepts in studying stationary stochastic processes. 
For our purposes here we need only consider these concepts 
as applied to discrete energy Signals. 
The autocorrelation function of any energy Signal b 


t 


is defined as: 


oo 


* 
pe tert 


* 
where the integer t is the time-shift and by is the compliex- 
conjugate of bi. If we only consider real-valued discrete 


energy Signals, the formula is simply: 


co 


bpp lt) a og rege : 


When t = O, we have 
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Da iCG)) = ve 
bb ‘ga 


which is equal to the energy E of the Signal, 


If we let s = t + T in the formula, we have 
Mop Ge 7 Peete tt i. cog 8 Sat 7 oe STS 7 ap bt D) 


Thus the autocorrelation function for a real-valued discrete 
energy Signal is an even function of the time-shift T. 

As a generalization of the autocorrelation function, 
the cross-correlation function of two real-valued discrete 


energy signals a, and by 1s defined as: 


foe) oO 


Pape = ert = Ny hh oc 


S=-& 


From the above formula it is easily seen that: 


> p6T) rc Ppal-T): 


Also we see that the autocorrelation function dbp Tt) is 
a special case of the cross-correlation UGG al eyi eel) 


occurring when the signals a, and b, are the same energy 


ie 


Siena]. 
6. Deconvolution 
DECouvelution is the inverse operation of convolution. 


Consider the problem that 


which means the discrete signal c, is equal to the discrete 


ie 


Signal a, convolved with discrete Signal b,- If we are given 


t 


Cy and a, and we wish to find b,, we need to deconvolve a 
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mrom CC, Lo get Die We can do the deconvolution both in 


t 
time-domain and in transform-domain, 
C 
Ci = Z ab 
Si 
C ay m tv 


The first few terms are 


lI 
rab) 
OY 
+ 
99 
oy 


CM 
tl 
99 
o 
oO 
+ 
9 
oi 
~ 
© 
oO 


Pietie first equation, we can solve Do» because we know 
both Cy and ao: After we get bo: Veomcanwuse bee seccong 


equation to solve b We can continue this iterative process 


1° 
until we find the signal b,. 

In the transform-domain, we know the convolution operation 
becomes a multiplication operation, Hence the deconvolution 
operation becomes a division operation. Thus if c = a *b, 
then after we take the Z-transform, we have 

CG) ee (2) B(Z ) 

If we are given Cy and ats we know C(Z) and A(Z). We can 

ie ahyaye! b, by first solving B(Z) which is equal to C(Z) divided 
by A(Z). Then we can transform B(Z) back to the time-domain 
to get the signal b,. 


To show the procedures stated above, let us do an 


example. Suppose 


yy 





Cc (Co, Cy» Co; Ca) = (12, 10, 14, 6) 


t 
a, = (a; a4) =a 2). 

We want to find b. such that 
cCc=a * b. 


Taking Z-transforms of Cc, and a,, we obtain 


12 + 10Z + 147z° + 6z° 


II 


C(Z) 
A(Z) = 4 + 22 


The Z-transform of Db, is 


2 3 : 
B(Z) = C(Z)/A(Z) = eee — eee ee at ees +- 37,7. 
Therefore 
Pe (Doubs, P,) = (3, 1, 3). 


One thing we must point out here is that when we convolve 
two finite-length discrete-time signals, the result is always 
a finite-length discrete-time signal. But when we deconvolve 
a finite-length signal from another finite-length signal, 
the results are not always of finite-length, For example, 


at 


C sae (Co; C1» Co» Cy) = (12, 10, 14, 8) 
a, = (ag, a,) = (4, 2) 
we will obtain 
by, = (by, by, be 1s) 
Co we ee 2 ee 41/8, =— 1/16, ...) 


Pen 1S Of infinite-length. This fact is related to the 
design of least error energy shaping filters, a subject we 


“It diseuss later on. 


18 





7. Stationary and Ergodic Stochastic Processes 


A stochastic process is a proces that assigns, 
according to a certain rule, a time function (continuous or 
discrete) f(t, ©) to every outcome ¢ of an experiment, It 
can be viewed as a function of two variables t and @. 

f(t, t) may represent: 
(1) A family of time functions (t and ¢ variables) 
(2) A single time function (t variable, ~t fixed) 
>) A random variable (t fixed, ¢ variable) 
(4) A single number (t fixed, t fixed). 

A stochastic process f(t) (we usually omit its 
dependence on ©) is stationary in the general sense if its 
expected value is a constant and its autocorrelation depends 
only on t,-to: 

Bent => Constant ( for any time t) 

E{f(t+t)f(t)} = R(t) (for any time t). 
mestochastic process f(t) is ergodic if its time averages 
equal its ensemble averages, 

We will assume the earth tilt data as an ergodic 
stationary time series in the design of single-channel 


predictor. 


cht, 





B,. THEORY OF SINGLE-CHANNEL PREDICTOR 
1, Method of Approach 
a. Inelegant Method 
By uSing the FFT (Fast Fourier Transform) we can 
determine the spectral density function of a time series, 
After we invert the spectral density function to the time 
domain we can get a Fourier series representation of the 
time series, We then extrapolate to future time to obtain 
a crude singe-channel prediction. 
b. Elegant Method 
The Wiener-Hopf equation can be used to design 
an optimum filter for a single channel predictor, The 
solution of the Wiener-Hopf equation is facilitated on a 
high-speed digiter computer by using the Levison recursive 
method which can be programmed. We will discuss this method 
of prediction in the following sections. 


2. Block Diagram of a Linear Time-Invariant Digital 
Filter 


The following is a block diagram of the input/output 


relationship of a linear-system: 








impulse 





+ response output: y 






sequence a 





x is the input signal to the linear-system, y is 
the output signal from the linear-system. The relation 


between x and y is: 
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where a is the impulse response sequence of the linear-systenm. 
The above equation is the time-domain relation, which 
Says the output y is equal to the convolution of the input 
x and impulse response sequence a of the system. 
In transform domain, the input-output relationship 


Simply becomes: 


Y(Z) = X(Z)A(Z). 
maac 1s, the Z-transftorm of the output signal y is equal to 
the product of the Z-transform of the input Signal x and the 
Z-transform of the impuise response sequence a of the system, 
A system that generates its output in this way is called a 
linear time-invariant digital operator, or a linear time- 
mevariant digital filter. The coefficients of the impulse 
response sequence are called the operator coefficients or 
tielLter coefficients. 

3. Least-Squares Shaping Filter 

From section 2.2, we know in a linear time-invariant 

system if we are given the input signal x and the impulse 
response sequence a of the system, we can get the output 
Signal y by convolving x and a. On the other hand, if we 
are given the input Signal x and we want the output signal 
y to be some desired form, how should we choose the impulse 
response sequence or filter coefficients a? Clearly this is 
a deconvolution problem. In transform-domain it becomes 


PeeciVvViSion problem. That is, 
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where A(Z) is the Z-transform of the desired filter coeffi- 
cients a. But we know from elementary algebra that even Y(Z) 
and X(Z) are both of finite terms, A(Z) may not be always 

of finite terms. In fact in most cases A(Z) will be a 
polynomial of infinitely many terms. This is to mean we 

want a filter with infinite-length. If in the real-world 
Situation we can only build a finite-length filter, then we 
can not obtain the exactly desired output Signal y. The 
difference between the actual output and the desired output 

is the error. For a specific finite-length filter, the optimum 
one is that one which will yield an output with smallet error. 
Like the regression analysis in statistics, here the smailest 
error means the least squared error, If we denote the desired 
SMepit Signal to be y and the actual output signal to be y, 


we want to minimize 


iiesoopcimum filter which can do this job is called a least- 
squares shaping filter since it can generate an output with 
least Squared error to the desired shape of output. 

4. Discrete-Time Wiener-Hopf HEquation 

In this section we will discuss how to determine 

the least-squares shaping filter's coefficients. Let the 
input Signal be: 

x = (Xo, By: x) 


let the desired output Signal be 


Mie Sipe ig 9 9 Vente 





let the filter coefficients be 
Cee, f_) 


The actual output signal y is the convolution of the input 


signal x and the filter coefficients f, that is: 


N 
veoe- x © f, 
Now we want to choose f which will minimize 


n+m Hat 
D) m 


A a 2 
ee Ve) >t Cy, - ff 
t=0 Cae Ss= 


The quantity I is minimized if its partial derivatives 


tiene respect to each of the coefficients fo: fi. eS ope aa 
equal zero. Therefore we want: 


Simplifying, we obtain: 


n+m 


mM 
eyo xe, + 2 Cpe he a) OX = 0 
t-o © *I t=0 s=0 § tS 


or 
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If we use the notation 


and 


where T 4s is the autocorrelation (for index j-s) of the 
input signal x and g. is the cross-correlation (for index j) 


of the output Signal y and the input signal x, we obtain 


After writing out this set of mtl equations, we have 


tg) 0 Sage) Se en a 20 
Por, + Liv, + for_y + + ae a ~ fy 
Pot t 47, + ToT * ei ee So 
for n + tiv +- for _o + + £0 = Ems 


Or in matrix form 
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O -—1 -~2 
a1 AO sai 
To ry ro 
rm amet m- 


= 
>) 


pod 


? 
= 
No 
1 


The above set of equations is called the normal equations 


and is the discrete-time Wiener-Hopf equation, 


5. A Computer Program for a Single-Channel Predictor 


In a least-Squares Shaping filter system, 


if the 


desired output signal y is equal to the input signal x 


shifted Ee time un 


Y, ~ % 
t ttt, 


its ahead, that is 


then the filter becomes a predictor since its objective is 


to predict the val 


the discrete-time 


m 
Es *s *5-s 
because 
n+m 
e 5 7 on 5: xt 


ue of x be 


Wiener-Hopf equation 


S (aly OS eae a HN, 

Hee (Gl = One iee , m) 
n+m 

.- = Fy xXx oe sale. 

-J uae Cet t-J alice 
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time units ahead, 


In this case, 





Writing out the equations in matrix form, and by noting 


a Ts for a real-value signal, we obtain 


O 1 2 m it 
O 
a Oo 1 El eet 
ma i 0  “m-2 = eee 
Um pmed man °°?’ © 76 | vt +m 





Although we must solve the above m+l1 equations to find the 


filter coefficients fo: tl. SS 


of the special form of the autocorrelation matrix 


m? We can take advantage 


ro a i, _ Um 

ar a a6 My Tn-1 
Sa | a6 Pm~2 

rn aL Pn-2 a 


in order to reduce the computational work. 
The autocorrelation matrix has the property that all 


the elements on any given diagonal are the same, It is also 
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a symmetric matrix. Hence the entire autocorrelation 
matrix R of order m+til does contain only mtl distinct 


Sements, namely =r oe ae Let us now see how 


0’ a 3 ; 
we can take advantage of this structure. 
The method we use 1s called the Levinson recursive 


method, Suppose that we want to solve the following matrix 


equation: 
| geen tg oe if fox 
Ol Zz rn 0 O 
cite, oa Tn-1 a By 
oy 27) is ee fo 7 Eo 
iene wal OC a em | 
where the autocorrelation coefficients To: Ta: To: alae en 


as well as the cross-correlation coefficients EQ: 84> Bo» 
"++,8, represent Dae mOWd GUuanbicies, ana the filter 
coefficients fy: i> Po: 


quantities. The recursive procedure solves the equation in 


a Ea represent the unknown 


a step-wise manner. The step k=0 is given as an initial 
condition, and then the steps k=1, k=2, ::-, k=m are done 
successively in a recurSive manner. The desired filter 
coefficients are the ones which result on the completion of 


the final step, namely step k=m. 


27 





Let us now show how we proceed from step k=n to step 
k=n+1. The completion of step k=n requires that we have 


computed the following quantities: 


and 


By definition, these quantities satisfy the matrix equations; 





6) uf 2 n+1 nO 
Sal Opel “n ani 
2 0 | eee et : = (1) 
7 ann 
n+1 ane na! “0 cs 
and 
9 Si Se aioe sl Fn0 Eo 
ry ro ry ees By 
“8 a 10 ees = 
: ii 
nese © n= 1 FO Yn 








Because of the special structure of the autocorrelation 
matrix, equation (1) may be written in the following 


equivalent form: 


“6 act re rntl u 8, 
= “OQ “al a: Ann 
To r, Yo Yue : | 
: (3) 
and O 
wee ne, n-l 86=6 TO 270 hy 


Let C. be a constant which we will determine in order to get 
a desired result. The desired result is after we multiply 
equation (3) by C, and add the result to equation (1), we 
want it to be identical to a new equation. That is, we want 


equation (4) and (5) to be identical. 


(6 pale Pnt+1 0) Es n n 
mi | Od en anit * Sn@nn O 
2a x0 “n=1 = 
| . (4) 
a foe a 0 
nn n nj 
ante ne nel "0 ce .) mn"n0 | Pn “nn | 
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0 gee nt+1 *n+1,0 n+1 
a Or 1 el anti, 1 0 
es i 0 ened 
‘ = : (5) 
Antiin ! 
r re © r Aa O 
+1 -1 0 n+1n+1 
n non Jn L | Q al 
In order for equations (4) and (5) to be identical, we 
require that 
On “ Cnn Ontl 
0 0 
; = , (6) 
0 0 
= aE Coo, | O | 
and 
2a ve an+1,0 
Ant 7 “nnn cel ol 
= CO), 
a + C4an1 ein 


or - on nd AntL ntl | 





From equation (6), we can determine the value of Ch by 


B +e a. = 0 
n nn 
or 
= By 
C me ee SS 
n a 
n 


After Ch has been determined, we can compute all the 


following quantities: 


ra a a a at a ¢ 
fei’ “nt1)O0? “n+1,1° “n+t1,2’ » Fn4d ntl 


+ Q@ + eee 


Baad = 2n4+1,0 nt1,0° n+2 Ale te ad if ant1,n+171 


then we get the counterpart of equation (1) of step k=n in 


step k=nt+1 as follows: 


0 1 2 n+2 4n+1,0 n+1 
ce Oo 1 °° 7 ntl Sree 
T9 ao 0 rn 4 ° 
hh ab 8) 
adel, AIpecal : 
ae? ntl ni m0 : Pn | 


Now, rewrite equation (5) as follows; 
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2 S58 Sel neki. wageal 3 
gee Oma Th nt+1,n : 
eo 7 Tn-1 
= (9) 
| net i 0 
eo OM en+1 0 | Feat 


bet us multiply equation (9) by a constant a as yet under- 
undetermined, and then add the result to equation (2). We 


obtain 


Pnt1 SOM onent1 ntl Et? 
rn Sal oa eet ae 
ei ; 
= (10) 
am ra onal ph aa 





We want equation (10) to be identical to the following 


new equation: 





T1 72 rie Fat .O 50 

FOr vn fol al By 

iL 36 eel ; _|- (11) 
—_ En 

i net Fo a | n+9., nt ent 





This means we require 


a Eo 
S45 Sy 
= GZS 
En Bn 
Be nad | Sn1 | 
and 
SAO aed. rel Le iO 
Ft nant Asal dl melt 
= ' C13 } 
nn or eal al Tatlin 


After qd, is substituted in equation (13), we can compute the 


following quantities: 


eco cui nie? 9? ned nei 
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Next, let 


Ler ol - f 


ec otf 


"pill QE mee Peoria eee petal ELS AL 


We obtain the counterpart of equation (2) of step kK=n in 


step k=ntl1 as follows: 


SOs |  n+2 ele 50 
5) cat eas I ened REL, al By 
rat 7 ¢ vn 
= (14) 
f g 
nel ned n+1 


vnt2 vntl aT a TO ‘i a 0 | | Yn 


We have completed the recursive procedure from step k=n to 


step k=nt+t1, To start at step k=0, we must determine 


and 


Foo’ Yo 


which satisfy the matrix equations: 


= (15) 


and 
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ry lo O 


4 


In equation (15), we set a 


we obtain 


Bo = Vy: 


From equation (16), we obtain 
g 
0 = 
00 2, 2 Me 00s 


Starting from step k=0 as above, we 


Meesteps k=1, 2, 3, --:, m, At the 


00 


. (16) 


Yo 


= 1 and solve Oy» B 


can proceed recursively 


final step k=m, the 


values obtained for the filter coefficients, namely 


will be the required solution. A computer program based on 


the above recursive method is included in Appendix A. 
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Pe eee un eRACLERISTICS OF FIDAL DATA 


A. INTRODUCTION 

By Newton's law of universal gravitation, each body in 
the universe attracts another body with a force directly 
proportional to the product of their masses and inversely 
proportional to the square of the distance between them, 


We may express this fact in a formula: 





1 2 


masses of the two bodies in the universe and r is the diStance 


where Ms is the gravitational force, m, and m, are the 


between them. Gis called the universal constant of grav- 


meat1On and 1S approximately equal to 6.67 x ‘LO none ee 
Celestial bodies such as the earth, the sun, and the 

moon are, of course, obeying the law of gravitation. Because 

they are not point masses and have large dimensions, different 

points on the surface of the body will be affected by differ- 

ent forces. Tidal forces are generated due to this effect. 

In this chapter we will give a general description of the 


emMaracteristics of the earth tidal data. 


B, GRAVITATIONAL POTENTIAL AND TIDAL POTENTIAL 
1. Basic Formulas From Mathematics 
Since in the derivation of tidal potential we will 


use some formulas from applied mathematics and spherical 
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trigonometry, 
section, 


a, 


Legendre Polynomials 


we give a brief discussion of them in this 


The Legendre Polynomials are defined by: 


2n- So). ss 
P(x) = (ee Lise 3 ) 1 


+> 


We observe that P(x) 1S a polynomial of degree n. 


few Legendre polynomials are 


Po(x) = 1 P(x) 
P,(x) = xX P,(x) 
Py(x) = 1/2(3x7-1) Pa (x) 


The Legendre polynomials can 


formula: 
1 





P(x) = 


After we have obtained Pox) 
Ziven above, we can also use 


formula to obtain the higher 


2nt+1 
n+1 





P41 6*) = xP (x) - 


a 


Nf1-2xt+t? 


for Legendre polynomials. 


The function is 


n(n-1)(n-2)(n-3 
2.4(2n-1)(2n-3 





= n-4 e 


as follows: 


3 


I 


1 /2(5x 


4 


IW oy Gero.< 


: 
1/8(63x° 


I 





Lo 


The first 
-3x ) 
-30x7+3) 
3} 
—-7Ox + 15x). 


also be expressed by Rodrique's 


1)", 


and P(x) by the formulas 


the following recurrence 


degree Legendre Polynomials: 


A 
n+1 


Bea) 


called tne generating function 


we can prove the following result: 
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By using the binomial theorem 





7 
a 


= al eee = y Pp (x)t™ ; 
n 
\| 2 2 
1-2xtt+t 
b, The Cosine Law for Spherical Triangles 
In a spherical triangle, if we denote the three 
sroaes tO be a, b, c and the three angles to be A, B, C; the 


following relations hold. 


Gos a = cos b cos c + sin b Sin c cos A, 


cos A = -—- cos B cos C + Sin B Sin C COS a. 


We call the above formulas the cosine law for spherical 
triangles. 


2. Der val lon GOltuGreavrbautlonal Potential and Tidal 
Potential. 


a. Gravitational Force, Gravitational Field 
Intensity, and Gravitational Potential 


According to Newton's law of universal gravita- 


mien, 2 point mass m, at a distance r from another body of 


O 


mass M will experience a force 


> MoM > 
F = -G —o i. 
5 te 
-> > 
where 1. 1S a4 unit vector in the direction of r. i is 


@alled the gravitational force. The minus sign denotes it 
imine abtractive force. 

lmstineombhemerayitabronal £1e¢ld intensity at a 
Heit to be the gravitational force experienced by a unit 


Masso atetnat point. if the unit mass is at a distance r 
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from the body with mass M, we will have 


where I denotes the gravitational field intensity. 
The gravitational potential - Ave pote whien 


lies at a distance r from a body with mass M is defined as 


follows: 
ee gd i Sat 
& ie 5 r r 
= J -G a dr = -G uw 
r 
a r 


Another way to show the relation between Ss and i is 


I = -VV_= -grad V 

g g g 
Since Ug is a scalar while Fe and Le are VeCerors - Me Ss 
much easier to handle than _ and I, In deriving the 


theoretical tide on earth induced by the mcon and the sun, 
we will use tidal potential instead of tidal force. 


b. Gravitational Potential and tidal potential at 
a Point on the Surface of Earth 
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Figure 1 


In Figure 1, let E be the earth, P is any point 
on the surface of E,. The gravitational potential at P created 


by a celestial body of mass M is 
= M 
V = -G R 


where R is the distance PB. Notice we are assuming that the 
whole mass of the celestial body is concentrated at a point 


B. By the cosine law for plane triangles, we have: 


R? = at F ee - Zgarcos €, 


Therefore we can rewrite the gravitational potential V as 
follows: 
V = -G u 


XG + r? - 2ar cos € 


=) -1/2 
= _¢ ! 1-28) cos € + 3)" 
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After we expand the above equation, we get a more compact 


form: 
E Pi(cos €)(=) 


where P_(cos —) are the Legendre polynomials, The first 


term of V is 


Wibeen 1S a constant and has no contribution to the gravitational 
force. Therefore it has no physical significance and we will 


disregard it, The second term is 


- «<M a. ig M 
V5 G = (cos &) G 2 ae COSuic., 


To evaluate the contribution of V, to the gravitational 
force, let us choose a cartesian frame of reference centered 
ere, the x axis running in the direction OB and the y axis 
running in the direction OO'. With this choice of coordinate 


system we have 


a cos € = x ; a Sin € = y : 
poe tne Contribution of Vy to the gravitational field intensity 
is 

‘= = _M = _M a 


>, ? : A . : : 
where 1, 1s a unit vector in the direction of xX axis, This 


shows that iy LS wUMeEOrMuctavyltacrtonal field. It has a 
constant magnitude G So and a constant direction which is 
ie 


parallel to the xX axis. Hence Vy will not create tidal effect. 
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mine third term of V is 





Weer 1/23 cos-e — 1) (2 : 
2 r O aa ) =) 
= -G aa: 1/j2 (3acos“é - a) 
1 i 

= -G "3 +1/2 (3x* - x” - y*) 
r 

= -G 4-1/2 (2x7 - y*) 
r 


in our choice of coordinate system. The gravitational field 


intensity created by Vo 1s 
> 
eae ye | =G —— (2x* - y*) 
2 2 3: 
PAD 
M a a 
= G $3 C2 te cy 7) 


> = 
where i, and Sy are unit vectors in the direction of x axis and 


+ 
fasts. We observe that I, is different at different locations 


2 
because it depends upon x and y. This non-uniform gravita- 
tional field intensity will create a tidal force at all points 
on the earth. Hence Vo is responsible for the tidal effect. 


We can also show that after V every term of V is responsible 


O°! 


for the tidal effect. Therefore we call 


ie encRreee Gco> =) (5) 


the tidal potential. 


eure Oventtalmol tmeuouriace Of HKarth Induced by 
the Moon and the Sun 


Although theoretically every celestial body in the 


universe creates a tidal potential on our earth, they are so 
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small that in practice we can neglect ail of them except 
those created by the moon and the sun. If we denote the tidal 
potential at a point P on the earth induced by the moon be 


Vw? by the sun be Vor we have 


Mn 00 at 
V.o=-G— z P (cos & ) (— ) 
m a n m ie 
m n=2 m 
Me 00 = Oe 
Vo. = -G — 2 P (eos €_) (—) 
S 1e i n S if 
S n=2 S 
where 
= " Wie 
Mn = the mass of the moon = 7.38 x 10 ‘kg, 
: 30 
Me Zaniomiacs cOtmt ne Sun — t,99l x 1Ow ke, 
a the distance from the center of the earth to the 
center of the moon, 
eae the distance from the center of the earth to the 
center of the sun, 
on = the zenith angle of the moon, 
aa = the zenith angle of the sun, and 
a = the radius of the earth = 63695 km. 


If we use the mean value of ae and the mean value of ee to 


calculate — and — , we obtain 
S 
ese 1<66 x 10“ 
ig O 


m 3°*84x10 


o 


i 
Il 


BEES = AOR ee 
1°49x10° 


43 





a a 
Since a and >— are both very small, we can express sis ana 








S 
Bie approximately as follows: 
Mn q 2 
V, 2 -G =" Py (cos &) (2) 
m m 
a 
= Gm . 1/2 (3c0s” 1 
See /2 (3cos™ & ~1), 
m 
Mo 7 2 
See ct Po (cos €_) (--) 
S S 
M_a® 2 
= -G ao ieeescOs &. —1). 
a S 
S 


For the same value of om and So? the ravionis 


= 2°17 


This indicates that the tidal potential due to the moon 
Zmserousity twice that due to the sun. Therefore the moon 
should be the preponderant source of the tidal effect on 


mie Carini, aS is indeed observed, 
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eee tHe PERTODICITIES IN THE TIDAL POTENTIAL INDUCED BY THE 
MOON AND THE SUN 


CELESTIAL POLE 


ZENITH 


f) 





Figure 2 
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In Figure 2, let P be a point on the earth with latitude 
L and longitude A, Q be the geographical point of the 
celestial body M which has declination d and local hour angle 
t at the instant shown in the figure. If the celestial body 
is the moon, then by the result of the last section we can 
write the tidal potential at P induced by it as: 


M a? 


a cell R 2 
ea = -G 3 Pieces COS &=—1)):, 


m 





Now we will derive the periodicities present in the tidal 
potential. By using the cosine law for spherical triangles, 


we can rewrite Vie. as 








2 
Ma 
VV = -G = . 1/2 [3(cos cos] 4a) Sin Gsin 7 cos t)?-a| 
5 r 
m 
Mae 5 
= -G 5 © 1/2 [3¢sin IL sin d + cos L cos d cos t -1| 
r 
m 


After we expand the above equation and algebraically 


Simplify, we obtain 





m eo [3(sin?1 mei oy (ein d= 13) 
r 
m 


<< 
It 
I 

Q 


+ sin 2L sin 2d cos t + Boece Poona cos 2t| 


K [3(sin®L-1/3) (sin?a-1/3] 


2 


+ k (Sin2k Sine dweeos. & ).. + K 6 Cos L cos“d COs aay) 


where 
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M a 


wi m ; 
K = -3/4 G 3 Loma COnsvant. 
Tm 





ters the latitude of P and is also a constant. d and t are 
the declination and local hour angle of the moon, respectively. 
They are changing from time to time. t will go through a 
complete cycle of 360° in one lunar day which is about 24 
hours and 50 minutes on the average. d will go through a 
complete cycle in about one month. Since in one lunar day 

d will change very little, we will consider the term 

K(sin 2L sin 2d cos t) as the tidal potential component with 
frequency 1 cycle/lunar day and the term K(cos“L cos“d cos2t ) 
With frequency 2 cycle/lunar day. We also call them diurnal 


~ 
and semi-diurnal component. The first term K|3(sin“L-1/3) 


(sin7a-1/3)| contains iim a which is an even function of d. 
Since d goes from zero to positive maximum and then back to 
zero and then goes to negative in a symmetric way, k|3(sin71 
-1/3) (sin?a-1/3) | will be a component with frequency i cycle 
mer Half month. We call it the fortnightly component. 

By the same procedure we can express the tidal potential 


at the point P induced by the sun as 


S 


V. = K. [3(sin*L-1/3)(sin?a-1/3) | 


a: Ki (sin2k Salim GeO sh bo). ot K_(cos“L cos“d cos2t ) 


M a“ 
S 


3 
r 
S 





where K. = -3/4 G 
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d and t of the sun go through a complete cycle with periods 
One year and one solar day, respectively. Hence the three 
terms of We have frequencies 1 cycle/half year, 1 cycle/solar 
day, and 2 cycle/solar day. 

If we do not neglect the higher order terms in the tidal 
potential, and also consider the quite complicated relative 
motion of the moon and the sun with respect to our earth 
in detail, the spectrum of the tidal potential should be 
composed of an infinite number of terms. But for practical 


purposes it is quite enough to consider only the most predom- 


inant nine terms, Table 1 iS a Summary of these terms. 
FREQUENCY 
NAME PERIOD (Chr) (degree/hr) FREQUENCY (rad./hr) 
K1 0.23934464E 02 0O.15041068E 02 O0.26251602E OO 
O1 OF ZeceLI330H O02 0.13943035E 02 0.24335176E OO 
Pl 0.24065887E 02 O.14958931E O02 0.26108247E OO 
QI. 0.26868347E 02 O.13398661E 02 0.23385066E O00 
M2 0.12420603E 02 0O.28984100E O02 0.50586778E OO 
S2 OFZ Oo0OCOE OZ O.30000000E 02 OF aZ255I98550H 00 
N2 0.12658348E 02 OP2S4397Z26H O02 0.49636668E O00 
K2 0.11967240E 02 OR se0e2 258.02 0.52503186E O00 
M4 OReZrosOl4y OL 0.57968201E C2 0.10117350E O1 


TABLE 1. Summary of Nine Principal Tidal Lines. 
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D. MISCELLANEOUS NOISE IN THE REAL TIDAL DATA 
Meteorological activities such as the rainfall, diurnal 
and seasonal temperature variations and pressure changes, 
timroawee NOlSe in the tidal tilt data. Other factors such 
as the ocean loading effect, daily and weekly traffic patterns 
also introduce some noise. However, to a lesser degree, 
these effects can be considered stationary and predictable 
over short periods of time. The effects due to local earth- 
quake activity or underground nuclear tests would introduce 
excessive noise into the real tidal tilt data. These noises 
are random in nature compared to other noise. We will 
remove these noises from the real tidal data when we use 


them as the input data of the single-channel predictor. 
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IV. APPLICATIONS OF SINGLE-CHANNEL PREDICTOR 
A, DISCUSSION OF THE PARAMETER'S VALUE USED IN THE COMPUTER _ 
PROGRAM 
Pepore using the "PREDICT" computer program, the first 
step is to convert a stationary continuous time signal to 
a discrete one by the method of Sampling. The sampling 


period T should be chosen according to the following formula: 


1 


2f 
m 





T < 


where fn is the highest frequency component of the time 
Signal. After we have obtained the discrete time signal, the 
next problem we will consider is how to choose the length 

of the input data; how far ahead should we predict and how 
long the filter length should be. Normally we should first 
do a spectral analyssi of the signal; find out the highest 
and the lowest frequency components. Then we should choose 
the length of the input signal data greater than the longest 
period of the spectral components. For the predicted 
distance, it is a good rule to choose it less than one- 
twelfth of the longest period of the spectral components and 
at the same time close to a small multiple of the shortest 
period of the spectral components. For the earth tidal tilt 
data, we have shown in chapter III that the most important 
spectral components are semi-diurnal, diurnal and fortnightly. 


Hence we choose the sampling period as one hour which is 


00 





considerably less than 1/2t_ (about 6 hours). Also, we 
choose the length of the input data to be 720 hours which 

is a little more than twice the period of the fortnightly 
components. The predicted distance we use is 24 hours which 
is less than one-twelfth of a fortnight and at the same time 
almost equal to two periods of the semi-diurnal components 
and one period of the diurnal components. 

Another parameter used in the "PREDICT" program is the 
length of the filter coefficients. It is the most important 
parameter which affects the computer CPU time used to do all 
the computations of the program. We will discuss it in 
detail in the next section. 

B, AN INVESTIGATION OF THE FILTER LENGTH USED IN THE 

"PREDICT" PROGRAM 

For a signal of pure sine wave form, it iS easy to show 
Ghat a iilter length of two is enough to predict the future 


data. From the basic trigonometry, we have 


Sain teeoes © tocos t sin T (1) 


Sin (t+T) 


Ssinee COS.) = cos t sin T (2) 


Sin (t-T) 


By adding equation (1) and (2) together, we obtain 


Hamer i ete sine cb—-!yo— 2 cos T sin t (3) 


or 


Secret = 2 cos LF sin t + (-1) sin (t-T) (4) 
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Equation (4) means we can use a two-length filter 


coefficients 
(2 cos T, - 1) 


to predict one unit sampling period T ahead if the signal is 
a pure sine wave. More specifically, let us suppose there 

is a sine wave signal with maximum amplitude 1 and sampling 
at a uniform time interval T which corresponds to 30 degrees. 


The discrete time signal will be 


i, 8, 2, 0, 


3 
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The two-length filter coefficients in this case are, 


(2 cos 30°, -1) = ( 3, -1). 


If we convolve the signal with the filter coefficients, we 


obtain a new signal 


1 V3 


which is just the original signal shifted one sampling period 
ahead. If we want to predict two sampling period ahead, we 
can obtain the required filter coefficients in a recursive 


way as follows: 


O24 





Sin..( t+2T ) Zee cme sa ne Ctl). + (—] ) sine t 


2 cos T [2 cos T sin t + (-1) sin (t-T)] 


fe (-1) Sdn t 


(4 cos“ T-1) sin t + (-2 cos T) sin (t-T) 


which in our example will be, 


aeos 30°-1. -2 cos 30°) = (2, - ¥3). 
Mmeeve COnvolve the original signal with (2, - 43), we will 
obtain a Signal, 

c 13 1, (3 ‘ 0 


ee y 9° 9” , 9 
4 


which is the original signal shifted 2 sampling periods 


N31 
-1, -"5 - 5, 0 


<a 


2 
2? ? 


d 


ahead. 

When a signal's spectrum contains more than one 
frequency, the sampling period T will correspond to different 
angle measurements for different frequency components. If 
we want to predict such signals one sampling period T ahead, 
the filter length in general will not be of finite length. 
With a finite length filter we can only approximately predict 
the time shifted signal. 

The earth tidal tilt data contains nine principal’ spectral 
components and various sources of noises. To predict the 
eretica! tiLit data by using the "PREDICT" program, the 
best way to choose the filter length is by doing a statis- 
tical analysis for some samples and then make a decision. 


Appendix B contains a summary of test results of RMS error 
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woe a iver Length used in the "PREDICT" program, Seven 
sample runs have been made. The conclusions are as follows: 

(1) The RMS error of the last 24 hours is much smaller 
than the RMS error of the whole 720 hours. 

Ozye, the filter length at 720, which is the maximum length 
for a 720 input data, causes larger prediction error 
than orner filter lengths. 

(3) The filter length between 60 and 660 does nat have too 
much eifeet on the predicted error, with only a little 
improvement in the central region. 

Based On the above results, we think that normally a filter 
length of 120 is enough for the prediction. But if we are 
Gimese concerned about the accuracy of the prediction, we may 
Meera filter length near 360, 

Meche length of input data is not 720, the rule of 
@mmoecing Lilter length will be 0.2 to 0.5 of the input data 
length. The length of the input data should be greater 
than the period of the lowest frequency component of the 


Signal. 


C. USING THE SINGLE-CHANNEL PREDICTOR AS AN INTERPOLATOR 
If for some reason there iS a missing data segment in 
an otherwise continuous record, we can apply the "PREDICT" 
mueocram to prediet the lost data. Since theoretically the 
program can be used to predict in both directions, it 1S a 
Zeoa rule to predict the lost data from both the forward 


Girection and the backward direction and then take the average 





value as the predicted data value. It is easy to show that 
the expected mean squared error of the average value is less 
than or at most equal to the expected mean Squared error 
caused by one direction prediction. Suppose at an arbitrary 
point in time the prediction error from the forward direction 
is Ej» from the backward direction is Eo» the averaged error 
is (ES + Ey) /2, we have, 





Mee) fe) oy fa * Be 7 PF iY 
ete 2 2 2 
4 
poe 2 
= | (E, - Ey)” > 0. 
or 
2 

fee 2 2 [1 Eo| 
me 22 i 2 


The above inequality means whenever we apply the single- 
channel predictor as an interpolator, we should always predict 
from both directions and use the averaged prediction value. 
Appendix C contains the results of six test runs on the real 
mivieetalt data to illustrate the application of singe- 


channel predictor. 
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APPENDIX B 


RESULTS OF SEVEN TEST RUNS OF PREDICTION 


UNDER DIFFERENT FILTER LENGTHS 


DATA USED 


Oe Instrument Number O07 


2. Instrument Location 
Berkeley California 
Crey otate 
Longitude: 122 Vea 
Latitude: 37°52' .6 N 


3. Date of Data 
Year 1970 


Year 1971 


61 





Os 


LENGTH 


60 
120 
180 
240 
300 
360 
420 
480 
040 
600 
660 
720 


MONTH 


0,405 
0.395 
0.386 
0.380 
0.379 
0.370 
0,367 
0.367 
0,365 
0.357 
0.354 
0.349 


oy 


TES LeRUN ENO. 1 


RNS_ ERROR 


LAST DAY 


0,139 
0.132 
0,106 
0,092 
0,089 
OF097 
0.114 
OL 11 
0.121 
O {Lom 
0.143 
0.109 





TEST RUN NO, 2 


FILTER RMS ERROR 
LENGTH MONTH LAST DAY 
60 0,289 0.229 
120 0.281 0.238 
180 0,276 0.255 
240 0.274 0.258 
300 0.273 0.257 
360 0.273 0.264 
420 0.272 0.261 
480 0.271 0.272 
540 0.270 0.275 
600 0.270 0.283 
660 0.270 0.292 
720 0.270 0.302 


63 





TEST RUN NO, 3 


FILTER RMS ERROR 
LENGTH MONTH LAST DAY 
60 0,304 0.144 
120 0.301 0.135 
180 Oe 229 0.135 
240 0.290 0.109 
300 0.284 Oat 
360 0.282 0.142 
420 0.279 0.133 
480 O221S Ona, 
040 0.278 0.128 
600 On2 77 0.130 
660 OF Zt O71 31 
720 0,277 0.134 


64 





TEST RUN NO, 4 


FILTER RMS ERROR 

LENGTH MONTH LAST DAY 
60 0.393 0.465 
120 0,407 0.474 
180 0.412 0.472 
240 0.413 0.449 
300 0.413 0.458 
360 0.414 0.463 
420 aul 0.464 
480 Ape 0.454 
540 0.421 0.464 
600 0.419 0.456 
660 0.422 0.444 
720 0.428 0.615 
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TEST RUN NO, 5 


FILTER RMS ERROR 
LENGTH MONTH LAST DAY 
60 0.194 OES 
120 Om ga 0.113 
180 Omg? O.117 
240 0.200 0.118 
300 0.203 Orme 9 
360 0.203 Orel? 
420 0.205 0.122 
480 0.203 0.120 
040 0.203 0.121 
600 0.204 0.122 
660 0.205 0.120 
720 0.207 0.144 


66 





TEST RUN NO, 6 


FILTER RMS ERROR 
LENGTH MONTH LAST DAY 
60 0.409 Oe LOS 
120 0,414 0.136 
180 Odd 1. OF 239 
240 0.410 0.134 
300 0.410 0,126 
360 0.410 0.142 
420 0.407 0.135 
480 0.408 Oh tee 
940 0.378 ORES 
600 OF (0 OF oO 
660 0.376 0.186 
720 Oma 5 0.271 





—" 
a a 


s ooh 6008 


TEST RUN NO, 7 


FILTER RMS ERROR 
LENGTH MONTH LAST DAY 
60 0.368 0.174 
120 0,397 0.171 
186 0.355 0,153 
240 0.3599 0.168 
300 0.350 0.178 
360 0.399 0.171 
420 0.367 0.165 
480 0.307 0.164 
040 0.369 0.177 
600 0,369 0.194 
660 0.371 OSS) Ie 
720 0.378 0.425 
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APPENDIX C 


RESULTS OF SIX TEST RUNS OF DATA INTERPOLATION 


BY APPLYING SINGLE-CHANNEL PREDICTOR 


DATA USED 
1, Instrument Number 07 
2. Instrument Location 
Berkeley California 
Casi, State 


Longitude: 122°14'.1 W 
Latitude: 37°52'.6 N 
3. Date of Data 


Apri 2201971 
Nort esOe, LOL 
May om LOL 
May oe LO’ 1 
May 24, 1971 
May sone 1971 


71 





ACTUAL 
DATA 
VALUE 


34. 
34, 
Saye 
35. 
go. 
39. 
35. 
BO. 
sO 
34. 
34. 
34. 
34. 
o4. 
34. 
34. 
34. 
34. 
34, 
34, 
34, 
34. 
34. 
34. 


990 
290 
100 
260 
400 
440 
380 
280 
100 
200 
700 
990 
450 
420 
470 
980 
730 
830 
850 
820 
720 
980 
450 
Tol 0) 


Aine 22. 


FORWARD 
PREDICTION 
VALUE ERROR 
34096 -0.046 
Sed 5 LEE! -0.134 
35.304 -0.204 
35.450 -0.190 
Somos -0.122 
35.940 -0.100 
35.502 -0.122 
35.444 -0.164 
35.364 -0,264 
35.250 -0.350 
Sho) o JOSH, -0.459 
S39), UC -0.552 
Dao OO -0.650 
35.158 =0:7386 
39.247 -O.777 
35.340 -Q. 760 
35.418 -0.688 
35.463 -0.633 
35.404 -0.554 
30.2¢/ -0.457 
39.141 -0.421 
347, 900 -0.410 
34.894 -0.444 
34.832 -0.452 


IS) Fae 
BACKWARD 
PREDICTION 
VALUE ERROR 
34.540 0.410 
34,631 Oy Sa, 
34.806 0.294 
epee & 0.267 
Ja as 0.247 
35.257 ORS 3 
35,300 0.080 
35,2595 OF0Z 5 
35.146 -0.046 
34.994 -0.094 
34.821 10), J Zaak 
34.674 -0,124 
34.580 -0.130 
34.565 -0.145 
Sa Oe -0.165 
Soe O9 -0,189 
Sao LO -0.186 
35.056 -0.226 
35.163 -C,313 
BSG Cee -0,401 
39.232 -0.512 
S3g dLeH@) -0.610 
35.146 -~0.696 
Shape LSM -O.797 


AVERAGE 
PREDICTION 

VALUE ERROR 

34.768 Orto2 
34.877 Oy LIL: 
39.095 0.045 
518) 4 PAIL OPUS 
Soot 0.063 
35.398 0.042 
35.401 -0.021 
35,349 -0,069 
35.255 -0,155 
99.122 -0.222 
er weclene -0.290 
34.888 -0.338 
34.840 -0.390 
34.861 -0.441 
34,941 -0,471 
35,054 -0.474 
35,167 -0.437 
35,259 -0,429 
S18) ate): -0.433 
Oma So -0.429 
39.186 -0,466 
35.090 -0.510 
oO 0 -0,570 
34.984 -0,604 


RMS ERROR FOR FORWARD PREDICTION = 0.464 
0.338 
RMS ERROR FOR AVERAGE PREDICTION = 0.396 


RMS ERROR FOR BACKWARD PREDICTION 


V2 





ACTUAL 
DATA 

VALUE 

33.360 
33.430 
33.500 
33.510 
33.470 
33.460 
33.450 
33.480 
33.530 
33.590 
33.660 
33.690 
33.680 
33.590 
33.420 
33.250 
ee 110 
32,990 
32.930 
32.920 
32.920 
32.960 
070 
e170 


js 0) oop Ol es] 0) 


FORWARD 
PREDICTION 
VALUE ERROR 
58) e410) -~0.560 
bo. o00 -0.470 
B18) 5 ae -0.390 
33.890 -0.380 
SO 90 -0.420 
B18) olen’ -~0.470 
34.010 -~0.560 
34.130 -~0,650 
34.210 -~0.680 
34.240 -~0.650 
34.180 -0.520 
34.080 -~0,.390 
33900 -0.280 
oon 170 -0.180 
33.540 -~0.120 
33,3590 -~0.100 
33.240 —-0.130 
Sieg hee, -0.190 
S514 178 -0,240 
33.220 -0.300 
33.330 -~0.410 
33.450 -0.490 
£36) ONY, -0.520 
33.680 -0.510 


RMS ERROR FOR FORWARD PREDICTION 
RMS ERROR FOR BACKWARD PREDICTION 
RMS ERROR FOR AVERAGE PREDICTION 


od 1 
BACKWARD 

PREDICTION 

VALUE ERROR 
32.860 O,500 
33,000 0.430 
CoO LOw 0,400 
33.120 0 OD) 
S03) L740) 0,350 
poe” 0.340 
33.080 ORS ho 
33.040 0.440 
Sie: . ONO) 6.520 
33.000 QO.590 
33.020 0.640 
33,030 0.660 
33.070 OIE ILG 
Sie) JLaEO 0.480 
33,130 O2290 
S15 EI) On240 
Bye) OES) 0.100 
32.870 8)5 1E210) 
32,760 OE 7O 
32.650 0.270 
32.610 Orolo 
O22, 900 OF 60 
32.620 0.450 
32.710 0.460 


GS 


AVERAGE 
PREDICTION 
VALUE ERROR 
Soo) -0.030 
33.450 -0.620 
33.495 0.005 
33.505 0.005 
33,9095 -0.035 
33,929 -0.065 
33.545 -~0.095 
33.985 -0.1095 
S36) 4 e118) -0.080 
$16) SAU) -0,030 
33.600 0.060 
33.990 OF LSo 
Sie, alte. Oe hess. 
33,440 0.150 
33.335 0.085 
Son ZOU 0,020 
33,125 -0,0195 
615g OPAS. -~0.035 
O14 a eh eh) -~0.035 
S29 oO -O.015 
S20 10 -0.050 
33.020 -0.060 
$36) 5 JUGS, -0.035 
$16) 5 ICIS) -0.025 
0.435 
0.421 
0.072 





ACTUAL 
DATA 

VALUE 

32,240 
32.270 
32.340 
32.450 
32,080 
32.660 
32.680 
32.650 
B22000 
32.410 
oe LOO 
32.020 
31.860 
31.760 
ou, 770 
31,850 
31.950 
32.030 
32.080 
32.120 
32.100 
32.060 
32.010 
ei, 960 


FORWARD 
PREDICTION 
VALUE 
~323 
foUZ 
EAS, 
.618 
. 700 
Bikes, 
Oro 
oO 
.684 
.485 
5 PRS 
5 Lae, 
.080 
.042 
-050 


WO, 


a0 
ts) (6 
. 341 
. 369 


ERROR 


.083 
Raves 
SSS, 
.168 
5 LAW 
CaS, 
5 ee! 
~153 
kot 
.O75 
5 (OS: 
eboO 
.221 
.282 
. 280 
Ba), 
Zeyh 
mot 
.261 
. 249 
.200 
.241 
pA 
~185 


May 8, 


Wk 
BACKWARD 
PREDICTION 
VALUE ERROR 
32,006 0.234 
32,0695 O.205 
So loo OFT 
32,239 On 218 
32,320 0.260 
32.412 0.248 
32.478 0.202 
32.470 0.180 
BVA OOS: O.217 
32.162 0.228 
32.019 Oza 78 
31.845 O.175 
Si OS? 0.273 
31.594 0.166 
Sil A Sree 0.196 
31,621 Om2z20 
31.685 0.265 
oil (80) 0.290 
31.878 0.202 
SEA eialdl 0.179 
31.945 0.155 
312933 0.124 
3 LAOMS 07022 
31.902 0.058 
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AVERAGE 
PREDICTION 
ERROR 


VALUE 


32.164 
32.228 
S25 59 
32.428 
32,910 
32.600 
32.6495 
32.030 
Bs, See 
32.333 
Slow 
32,000 
31.883 
31.818 
31.815 
31.8695 
31.943 
32.033 
BIA, MONS) 
32.155 
32.147 
eylhy, Mat Fy 
32,069 
32.023 


RMS ERROR FOR FORWARD PREDICTION = 0.198 
RMS ERROR FOR BACKWARD PREDICTION 
RMS ERROR FOR AVERAGE PREDICTION = 0,045 


Oe Ales) 


S) (@ (SE 1S) (ee) 2) ei 2s, 72, 


0.076 
.042 
~O01 
.O21 
7070 
.060 
»035 
Ones 
.042 
.076 
RUSH, 
.020 
.624 
.O8 

.045 
.O15 
.007 
=O: 
.030 
.O035 
.048 
.057 
1OS0 
.063 





May 15, 1971 


ACTUAL FORWARD BACKWARD AVERAGE 

DATA PREDICTION PREDICTION PReprCrion 
VALUE VALUE ERROR VALUE ERROR VALUE ERROR 
Sa. 200 Je 30 -0.150 31.440 -0.240 SUL, Spe kD) -0.195 
ea, 200 oly Srae -0.150 31.570 -0.370 31.460 -0.260 
31.200 31.360 -0.160 31.670 -0.470 Sa UEsay ite. -0.315 
eri OO soo -0.200 opUR a8 -0.530 31,555 -0.365 
eal. i 70 31.450 -0,280 31.720 -0,990 SU ayers, -0.415 
o1.170 31.530 -0.360 31,670 -0.500 31.600 -0.430 
31.180 31.600 -0,420 ele @) -0.410 SEE he) s. -0,415 
51.230 31.660 -0.430 31.590 -0.320 31.605 -0.375 
31.3006 31.740 -0.440 31.550 -0.290 31.645 -0,345 
31.360 od, 7720 -0.410 31.570 -0,210 ol Ow. -0,316 
31.420 ill (AR, -0.290 Sil yt epe@) -0,210 31.670 -0,250 
31,430 31.570 -0.140 Sik ois l) -0.260 31.630 -0.200 
31,380 Sta 0 -0.030 BE Fe) -0.340 31.565 -0.185 
ou, 250 Mo 42119, 0.050 Si -0.450 31.480 -0.200 
31.100 31.040 0.060 31.6380 -0.980 31,360 -0,260 
2O2920 207590 0.030 31.5950 -0.630 Sis 2748) -0.300 
30.740 30.760 -0,020 31.370 -0.630 31.0695 -0.329 
30". 590 30.670 -0.080 )dk 6 7A a8 -~0.640 30.9090 -0.360 
30.500 30.680 -0.180 Sil, Jee -0.630 30.905 -0.405 
30.490 30.740 -0.250 31.020 -0.530 30.880 -~0.390 
30.520 30.840 -0.320 30.940 -0,420 B07 O50 -0.370 
30.640 30.980 -~0.340 30,940 -0.300 202700 -0,320 
30.810 3 ealeO -0.300 31.030 -0.220 242 O40 -0.260 
30.950 iL, AIRe, -0.260 31.180 -0.230 Le Lk She) -0.2495 


RMS ERROR FOR FORWARD PREDICTION = 0.261 
RMS ERROR FOR BACKWARD PREDICTION = 0.440 
RMS ERROR FOR AVERAGE PREDICTION = 0.321 
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May 24, 1971 


ACTUAL FORWARD BACKWARD AVERAGE 

DATA PREDICTION PREDICTION PREDICTION 
VALUE VALUE ERROR VALUE ERROR VALUE ERROR 
30,940 50.991 -~0.051 SE OLS 70.073 ol OO2 -~0,062 
30.930 oly Oke -0,089 Sn Ole, -0.082 Sly Ol -0.085 
30.980 31,112 -0.132 31.028 -0.048 3 O70 -0,090 
31.060 31.243 -0.183 31.067 -0.,007 31,155 =0..099 
31.180 SU Se) -~0.183 31.154 0.026 31.206 -0.078 
si. 300 31.446 -0,146 31.265 OF Ose 31.3595 -0.055 
ai 3 1 0 31.451 -0.081 31.367 0.003 S149 -0.039 
Siac Sh Le open | -0.011 31.452 -0.072 31.421 -~-0.041 
31.340 31.274 0.066 sia 30 -0,.090 31.392 -0.012 
er. 230 31050 0.140 el. Slee, -0.0895 31.202 (Oy Ox: 
60290 S0OnS9F OF Uo 31,140 -0.150 Si Ous -0.028 
30.750 30.671 0.079 20,902 -0,152 S10 Loo -0.036 
30,560 30,483 ONO 7a sO nOOg -0.109 30.576 -0.016 
507410 SO oo Oe Oils. 30.465 -0.055 30.428 -~-0.018 
30.330 30.359 -0.025 30,316 0.014 30,339 -0.0095 


SOs 1.0 30.358 -0.048 30RZo" ORCS, 30.304 0.0095 
30,340 S0n220 -0,080 30.235 O,1095 30.327 0.012 
30.460 30.544 -~0,084 30.294 0,166 30.419 0,041 
30.640 30.682 -0.042 30.416 0,224 30.549 O77 071 
30.780 30.804 -0.024 30.570 Oe 4A, 307 Ge7 0; O023 
30.880 S06 9 1 -~0,011 Z2UmbS & OF kos SIC) (el, Omg 
30.950 307,912 0.038 30.780 OrALee 30.846 0.104 
307090 30.904 OLGS6 302633 On low 30.868 Oral 
31,000 30.874 0,126 30.840 0-260 30.857 0,143 


RMS ERROR FOR FORWARD PREDICTION = 0.094 
RMS ERROR FOR BACKWARD PREDICTION = 0.121 
RMS ERROR FOR AVERAGE PREDICTION = 0,070 
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ACTUAL 
DATA 

VALUE 

30.300 
30.380 
30.440 
30,470 
30.470 
30.450 
30,420 
30.420 
30.440 
30.460 
30.930 
30.600 
30.630 
30.620 
30.540 
30.380 
30,210 
30.050 
ao. 930 
Zio. 900 
Zo.920 
ie. 990 
30.080 
30.180 


May 30, 1971 
FORWARD BACKWARD 
PREDICTION PREDICTION 

VALUE ERROR VALUE ERROR 
SCne 00m O0F090 30.380 -0.080 
30.380 0.000 30.450 -0.070 
30,340 0.100 30.490 -0.050 
30.320 0.150 30,510 -0.040 
30.320 0.150 30.510 ~-0.040 
30.340 0.110 30.460 -0.010 
30.390 0.030 30.350 0,070 
30.470 -0.050 30.250 6.170 
30.550 -0.110 30.200 DBA 
30.610 -0.150 30.190 one 70 
SOCOM OF lon | 30,220 0.310 
SORGCOm 0/0200" 30.260 0.3406 
30.530 0.100 30.330 0.300 
30.380 0.240 30.380 0.240 
30.200 0.340 30.420 ODO 
30,040 0.340 30.450 -0.070 
29.910 Ose OmeesO- 420 0 8©6- 0. 2110 
29.860 0.190 30.340 -0.290 
29.840 0.090 30.260 -0.330 
29.900 CHOOOMmESOMI ZO  -0.270 
SCR OLOME=OMOo Om aOrOo0| -0.170 
SO pdOm Onl lO GomOsom -0.090 
30.190 -0.110 30.090 -0,.010 
30.260 -0.080 30.140 0.040 


RMS ERROR FOR FORWARD PREDICTION 


RMS ERROR FOR BACKWARD PREDICTION 


RMS ERROR FOR AVERAGE PREDICTION 


U7 


AVERAGE 

PREDICTION 
VALUE ERROR 
S10 relore -0.085 
30.415 -0.035 
30.415 0.025 
30.4195 0.055 
30.415 G.Q59 
30.400 0.050 
30.370 0.050 
30,360 0.060 
30.375 0.065 
30.400 0.060 
30.430 O.100 
30.440 0.160 
30.430 0.200 
30.380 0.240 
30,310 0,230 
30.245 0.135 
30.1695 0.045 
30.100 -0,050 
30.050 -0.120 
30.035 -0.1395 
30.050 -0.130 
30.095 -0.105 
30.140 -0.060 
30,200 -0.020 
0,158 

- 0.199 

OFFS Bi Bes: 
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